6  Shape functions and finite elements

6.1 Shape functions

Linear shape functions

N = [ValueAtLF(-1), ValueAtLF(1)]
P = monomials(0:1, -1 .. 1)
M = [n(p) for p  P, n  N]
L11 = inv(M) * P

fplot(L11...)

Hermite polynomials

h = 2.0
N = [ValueAtLF(-1), DerivativeAtLF(-1), ValueAtLF(1), DerivativeAtLF(1)]
P = monomials(0:length(N)-1, IHat)
M = [n(p) for p in P, n in N]
H1 = inv(M) * P
fplot(H1...)

Hermite like polynomials with middle node

h = 2.0
N = [ValueAtLF.(0:h/2:h); DerivativeAtLF.([0, h])] |> vcat
P = monomials(0:length(N)-1, 0 .. h)
M = [n(p) for p in P, n in N]
H = inv(M) * P
fplot(H...)

Bilinear shape functions

WGLMakie.activate!()
V = [-1 1 1 -1; -1 -1 1 1]
N = [ValueAtLF(p) for p in eachcol(V)]
P = mmonomials(2, 1, QHat)
M = [n(p) for p in P, n in N]
L21 = inv(M) * P
fplot3d(L21)

Serendipity shape functions

V = [-1 0 1 1 1 0 -1 -1; -1 -1 -1 0 1 1 1 0]
P = mmonomials(2, 2, QHat, (p1, p2) -> p1 + p2 < 4)
N = [ValueAtLF(p) for p in eachcol(V)]
M = [n(p) for p in P, n in N]
Q8 = inv(M) * P
fplot3d(Q8)

Conforming rectangular plate element shape functions

V = [-1 1 1 -1; -1 -1 1 1]
P = mmonomials(2, 3, QHat)
N = vcat(
    [
        [
            ValueAtLF(p),
            PDerivativeAtLF(p, [1, 0]),
            PDerivativeAtLF(p, [0, 1]),
            PDerivativeAtLF(p, [1, 1])
        ]
        for p in eachcol(V)
    ]...
)
M = [n(p) for p in P, n in N]
H4 = inv(M) * P

fplot3d(H4)

6.2 Finite element stiffness matrices

@variables l, a, b, EA, EI, ν;

6.2.1 Truss

D(u) = 2 / l * u'
ae(u, δu) = EA * l / 2 * integrate(D(u) * D(δu), IHat)
Ke = simplify.([ae(s1, s2) for s1  L11, s2  L11])

\[ \begin{equation} \left[ \begin{array}{cc} \frac{EA}{l} & \frac{ - EA}{l} \\ \frac{ - EA}{l} & \frac{EA}{l} \\ \end{array} \right] \end{equation} \]

6.2.2 Beam

H1l = [H1[1], l / 2 * H1[2], H1[3], l / 2 * H1[4]]
D(w) = 4 / l^2 * w''
ae(w, δw) = EI * l / 2 * integrate(D(w) * D(δw), IHat)
Ke = simplify.(l^3 / EI * [ae(s1, s2) for s1  H1l, s2  H1l])

\[ \begin{equation} \left[ \begin{array}{cccc} 12 & 6 l & -12 & 6 l \\ 6 l & 4 l^{2} & - 6 l & 2 l^{2} \\ -12 & - 6 l & 12 & - 6 l \\ 6 l & 2 l^{2} & - 6 l & 4 l^{2} \\ \end{array} \right] \end{equation} \]

6.2.3 2D Poisson equation

Integration on the reference quadrilateral

D(u) = [2 / a * ∂x(u), 2 / b * ∂y(u)]
ae(u, v) = a * b / 4 * integrate(D(u)  D(v), QHat)
Ke = simplify.(expand.(6 * a * b * [ae(n1, n2) for n1  L21, n2  L21]))

\[ \begin{equation} \left[ \begin{array}{cccc} 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} & - a^{2} - b^{2} & - 2 a^{2} + b^{2} \\ a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) & - 2 a^{2} + b^{2} & - a^{2} - b^{2} \\ - a^{2} - b^{2} & - 2 a^{2} + b^{2} & 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} \\ - 2 a^{2} + b^{2} & - a^{2} - b^{2} & a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) \\ \end{array} \right] \end{equation} \]

Alternative derivation in physical coordinates

V = [0 a a 0; 0 0 b b]
N = [ValueAtLF(p) for p in eachcol(V)]
P = mmonomials(2, 1, QHat)
M = [n(p) for p in P, n in N]
L21p = inv(M) * P

ae(u, v) = integrate((u)  (v), (0 .. a) × (0 .. b))
Ke = simplify.(expand.(6 * a * b * [ae(n1, n2) for n1  L21p, n2  L21p]))

\[ \begin{equation} \left[ \begin{array}{cccc} 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} & - a^{2} - b^{2} & - 2 a^{2} + b^{2} \\ a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) & - 2 a^{2} + b^{2} & - a^{2} - b^{2} \\ - a^{2} - b^{2} & - 2 a^{2} + b^{2} & 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} \\ - 2 a^{2} + b^{2} & - a^{2} - b^{2} & a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) \\ \end{array} \right] \end{equation} \]

Nicely formatted output. TODO: Handle LaTeX-String from latexify

\[ \mathbf{K}^\mathrm{e} = \frac{1}{6ab} \cdot L"\left[ \begin{array}{cccc} 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} & - a^{2} - b^{2} & - 2 a^{2} + b^{2} \\ a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) & - 2 a^{2} + b^{2} & - a^{2} - b^{2} \\ - a^{2} - b^{2} & - 2 a^{2} + b^{2} & 2 \left( a^{2} + b^{2} \right) & a^{2} - 2 b^{2} \\ - 2 a^{2} + b^{2} & - a^{2} - b^{2} & a^{2} - 2 b^{2} & 2 \left( a^{2} + b^{2} \right) \\ \end{array} \right]" \]